# load required libraries
library(tidyverse)
library(langcog) # source: https://github.com/langcog/langcog-package
library(psych)
library(lme4)
# set theme for ggplots
theme_set(theme_bw())
d_all <- d_us_ad_pilot_22 %>% rownames_to_column("subid") %>%
full_join(d_gh_ch_22 %>% rownames_to_column("subid")) %>%
full_join(d_th_ad_22 %>% rownames_to_column("subid")) %>%
full_join(d_th_ch_22 %>% rownames_to_column("subid")) %>%
full_join(d_vt_ad_22 %>% rownames_to_column("subid")) %>%
full_join(d_vt_ch_22 %>% rownames_to_column("subid")) %>%
column_to_rownames("subid")
Joining, by = c("subid", "get angry", "choose what to do", "figure out how to do things", "feel guilty", "feel happy", "hear things", "get hungry", "get hurt feelings", "feel love", "feel pain", "pray", "feel proud", "remember things", "feel sad", "feel scared", "sense when things are far away", "sense temperatures", "feel shy", "feel sick, like when you feel like you might vomit", "smell things", "think about things", "feel tired")
Joining, by = c("subid", "get angry", "choose what to do", "figure out how to do things", "feel guilty", "feel happy", "hear things", "get hungry", "get hurt feelings", "feel love", "feel pain", "pray", "feel proud", "remember things", "feel sad", "feel scared", "sense when things are far away", "sense temperatures", "feel shy", "feel sick, like when you feel like you might vomit", "smell things", "think about things", "feel tired")
Joining, by = c("subid", "get angry", "choose what to do", "figure out how to do things", "feel guilty", "feel happy", "hear things", "get hungry", "get hurt feelings", "feel love", "feel pain", "pray", "feel proud", "remember things", "feel sad", "feel scared", "sense when things are far away", "sense temperatures", "feel shy", "feel sick, like when you feel like you might vomit", "smell things", "think about things", "feel tired")
Joining, by = c("subid", "get angry", "choose what to do", "figure out how to do things", "feel guilty", "feel happy", "hear things", "get hungry", "get hurt feelings", "feel love", "feel pain", "pray", "feel proud", "remember things", "feel sad", "feel scared", "sense when things are far away", "sense temperatures", "feel shy", "feel sick, like when you feel like you might vomit", "smell things", "think about things", "feel tired")
Joining, by = c("subid", "get angry", "choose what to do", "figure out how to do things", "feel guilty", "feel happy", "hear things", "get hungry", "get hurt feelings", "feel love", "feel pain", "pray", "feel proud", "remember things", "feel sad", "feel scared", "sense when things are far away", "sense temperatures", "feel shy", "feel sick, like when you feel like you might vomit", "smell things", "think about things", "feel tired")
Pooling all participants from all sites together into a common factor structure.
reten_all_PA <- fa.parallel(d_all, plot = F); reten_all_PA
Parallel analysis suggests that the number of factors = 4 and the number of components = 3
Call: fa.parallel(x = d_all, plot = F)
Parallel analysis suggests that the number of factors = 4 and the number of components = 3
Eigen Values of
Original factors Simulated data Original components simulated data
1 8.53 0.33 9.13 1.29
2 1.28 0.25 1.92 1.24
3 0.62 0.21 1.26 1.20
4 0.20 0.18 0.83 1.17
reten_all_par <- reten_all_PA$nfact
efa_all_par <- fa(d_all, nfactors = reten_all_par, rotate = "oblimin",
scores = "tenBerge", impute = "median")
convergence not obtained in GPFoblq. 1000 iterations used.
heatmap_fun(efa_all_par) +
labs(title = "Parallel Analysis")
Joining, by = "capacity"
Joining, by = "factor"
scoresplot_fun(efa_all_par, target = "all") +
labs(title = "Parallel Analysis")
Ignoring unknown aesthetics: y
scoresplot_fun(efa_all_par, target = c("ghosts", "God", "children")) +
labs(title = "Parallel Analysis")
the condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedIgnoring unknown aesthetics: y
itemsplot_fun(efa_all_par, target = c("ghosts", "God", "children")) +
labs(title = "Parallel Analysis")
the condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedJoining, by = "capacity"
|============================================ | 73% ~1 s remaining
|============================================= | 75% ~1 s remaining
|============================================== | 76% ~1 s remaining
|=============================================== | 77% ~1 s remaining
|================================================ | 79% ~1 s remaining
|================================================= | 81% ~1 s remaining
|================================================== | 83% ~0 s remaining
|=================================================== | 85% ~0 s remaining
|==================================================== | 86% ~0 s remaining
|===================================================== | 88% ~0 s remaining
|====================================================== | 89% ~0 s remaining
|======================================================= | 92% ~0 s remaining
|======================================================== | 93% ~0 s remaining
|========================================================= | 94% ~0 s remaining
|========================================================== | 95% ~0 s remaining
|=========================================================== | 97% ~0 s remaining
|============================================================ | 99% ~0 s remaining
Joining, by = c("factor", "capacity", "order")
reten_all_BIC <- VSS(d_all, plot = F); reten_all_BIC
Very Simple Structure
Call: vss(x = x, n = n, rotate = rotate, diagonal = diagonal, fm = fm,
n.obs = n.obs, plot = plot, title = title, use = use, cor = cor)
VSS complexity 1 achieves a maximimum of 0.88 with 1 factors
VSS complexity 2 achieves a maximimum of 0.92 with 2 factors
The Velicer MAP achieves a minimum of 0.01 with 3 factors
BIC achieves a minimum of -645.05 with 3 factors
Sample Size adjusted BIC achieves a minimum of -196.92 with 5 factors
Statistics by number of factors
vss1 vss2 map dof chisq prob sqresid fit RMSEA BIC SABIC complex eChisq
1 0.88 0.00 0.018 209 1964 1.7e-282 11.0 0.88 0.096 539 1203 1.0 2336
2 0.56 0.92 0.011 188 899 4.0e-93 7.5 0.92 0.065 -383 214 1.5 696
3 0.43 0.82 0.011 168 501 9.0e-35 6.2 0.93 0.047 -645 -112 1.8 265
4 0.43 0.81 0.014 149 383 1.5e-22 5.7 0.94 0.042 -633 -160 1.9 184
5 0.40 0.73 0.018 131 280 6.9e-13 5.4 0.94 0.036 -613 -197 2.1 134
6 0.40 0.73 0.023 114 227 1.6e-09 5.1 0.95 0.033 -550 -188 2.2 103
7 0.40 0.68 0.027 98 182 5.7e-07 4.8 0.95 0.031 -486 -175 2.3 75
8 0.37 0.60 0.033 83 134 3.3e-04 4.5 0.95 0.026 -432 -168 2.5 52
SRMR eCRMS eBIC
1 0.074 0.078 911
2 0.041 0.045 -586
3 0.025 0.029 -880
4 0.021 0.026 -832
5 0.018 0.024 -759
6 0.016 0.022 -674
7 0.013 0.020 -594
8 0.011 0.019 -514
reten_all_bic <- data.frame(reten_all_BIC$vss.stats %>% rownames_to_column("nfact") %>% top_n(-1, BIC))$nfact %>% as.numeric()
efa_all_bic <- fa(d_all, nfactors = reten_all_bic, rotate = "oblimin",
scores = "tenBerge", impute = "median")
heatmap_fun(efa_all_bic) +
labs(title = "Minimizing BIC")
Joining, by = "capacity"
Joining, by = "factor"
# labs(title = "Figure 1: Shared conceptual structure")
scoresplot_fun(efa_all_bic, target = "all") + #, highlight = "supernatural") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1,
color = c(rep("black", 8),
rep("#984ea3", 2)),
face = c(rep("plain", 8),
rep("bold", 2)))) +
scale_fill_manual(values = c("#e41a1c", "#4daf4a", "#377eb8")) +
labs(title = "Minimizing BIC")
Ignoring unknown aesthetics: yScale for 'fill' is already present. Adding another scale for 'fill', which will
replace the existing scale.
# labs(title = "Figure 2: Factor scores")
scoresplot_fun(efa_all_bic, target = c("ghosts", "God", "children")) +
scale_fill_manual(values = c("#e41a1c", "#4daf4a", "#377eb8")) +
labs(title = "Minimizing BIC")
the condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedIgnoring unknown aesthetics: yScale for 'fill' is already present. Adding another scale for 'fill', which will
replace the existing scale.
itemsplot_fun(efa_all_bic, target = c("ghosts", "God", "children")) +
labs(title = "Minimizing BIC")
the condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedthe condition has length > 1 and only the first element will be usedJoining, by = "capacity"
|======================================== | 67% ~1 s remaining
|========================================= | 68% ~1 s remaining
|========================================== | 70% ~1 s remaining
|=========================================== | 72% ~1 s remaining
|============================================ | 73% ~1 s remaining
|============================================= | 75% ~1 s remaining
|============================================== | 76% ~1 s remaining
|=============================================== | 78% ~1 s remaining
|================================================ | 80% ~1 s remaining
|================================================= | 81% ~1 s remaining
|================================================== | 83% ~1 s remaining
|================================================== | 84% ~1 s remaining
|=================================================== | 85% ~0 s remaining
|==================================================== | 87% ~0 s remaining
|===================================================== | 88% ~0 s remaining
|====================================================== | 90% ~0 s remaining
|======================================================== | 92% ~0 s remaining
|========================================================= | 94% ~0 s remaining
|========================================================= | 94% ~0 s remaining
|========================================================== | 96% ~0 s remaining
|=========================================================== | 97% ~0 s remaining
|============================================================ | 99% ~0 s remaining
Joining, by = c("factor", "capacity", "order")